Vortex viscosity in the moderately clean limit of layered superconductors 
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We present a microscopic calculation of the energy dissipation in the core of a vortex moving in 
a two-dimensional or layered superconductor in the moderately clean regime. In this regime, the 
quasiclassical Bardeen-Stephen result remains valid in spite of the strong correlations between the 
energy levels. We find that the quasiclassical expression applies both in the limit of fast vortex 
motion (with transitions between smeared levels) and in the limit of slow vortex motion (with 
nearly adiabatic dynamics). This finding can be related to the similar result known for the unitary 
random-matrix model. 
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I. INTRODUCTION 



At low temperature the energy dissipation in the vor- 
tex cores is the main source oLresistivity in the mixed 
state of type-II superconductorsu. If a supercurrent flows 
through a superconductor, it exerts a force on the vor- 
tices. Unless pinned by impurities or inhomogeneities, 
the vortices are brought into motion, which in turn leads 
to dissipation. Depending on the level of disorder, the 
vortices may move at different angles with respect to 
the direction of the supercurrent. At weak disorder, the 
vortices move together with the supercurrent ("ballistic 
limit", with the Hall angle close to ir/2). At sufficiently 
strong disorder, the vortex motion is directed perpendic- 
ular to the supercurrent ( "dissipativc limit" , small Hall 
angle). Both limits ace. jssjell jinderstood within the qua- 
siclassical descriptionalj'CrtMil A simplified approach de- 
scribing the dissipative limit goes back to the theory of 
Bardeen and Stephen treating the vortex as a region of 
normal phase inside a superconductorta. In spite of ne- 
glecting the structure of the quasiparticle excitations in 
the vortex coreEl, the Bardeen-Stephen theory gives the 
same result (up to a num<j»jka|Ufactor) as the accurate 
quasiclassical calculationaBalaQLl. 

It has been recently suggested that the microscopic 
structure of the core excitations may play a much 
more prominent role in that part of the dissipative 
regime where the excitation spectrum remains dis- 
crete (sharp quasipapticle levels), specifically in layered 
superconductorslljOllallj. In the clean limit (i.e., for 
scattering rates H/t much smaller than the superconduct- 
ing gap A), the motion of quasiparticles in the vortex core 
is ballistic: they cross the vortex core many times be- 
fore scattering off impurities. Therefore in this limit the 
spectral properties are sensitive to the details of disorder 
realization. In the superclean regime (H/t <§; A 2 /Ep, 
where Ep is the Fermi energy) the levels inside a two- 
dimensional TOftex split into two sets ('combs') of equally 
spaced levelsEil. The transformation of this correlated 
spectrum into a featureless uncorrelated one with in- 
creasing disorder proceeds in distinct steps: within the 
moderately clean regime (A 2 /Ep <§; H/t <C A) a new 



intermediate .region (A 2 /E F < H/t < A^A/E F ) has 
been foundtJO where the comb structure remains pre- 
served but is randomly shifted in energy (the number 
of impurities in the core has to be small enough to pre- 
serve the 'combs', while being large enough to random- 
ize the overall shift). Increasing the impurity concentra- 
tion further (already in the moderately clean limit _£pr 
weak impuritiestj, and ultimately in the dirty limitlij), 
the 'comb' structure is destroyecLwith a crossover to the 
class C random-matrix ensembleE£l. 

Keeping within the moderately clean regime with a 
discrete spectrum, in two-dimensional (or layered) super- 
conductors there further exist two limits of dissipationt^: 
The first one applies to the slowly moving vortex with 
a discrete quasiparticle spectrum ("discrete-spectrum 
regime") where the dissipation is due-jtjO,Landau-Zener 
transitions between individual levelsOiia. The second 
limit is with levels smeared into a continuum either by 
the vortex motion or by inelastic processes ( "continuum- 
spectrum regime"); the dissipation then is given by the 
linear-response Kubo formulalij. At low temperatures, 
the inelastic smearing can be neglected and the crossover 
between these two regimes is controlled by the vortex 
velocity with the characteristic velocity given by vk ~ 
{A/pf)/ \fkp~l (here, A is the superconducting gap and I 
denotes the elastic mean free patf£9). 

In the framework of random-matrix models with time- 
dependent Hamiltonians, the dissipation in the discrete- 
spectrum and continuum-spectrum regimes was consid- 
ered by WilkinsorJla. He finds that, in the case of the 
unitary Wigner-Dyson ensemble, the linear dissipative 
response remains valid in the whole range of velocities, 
both in the continuum-spectrum (high velocity) and in 
the discrete-spectrum (low velocity) regimes. Based on 
this fact and on the similarity between the unitary and 
class C ensembles, it was shown that with the class 
C level statistics the dissipation rate nearly follows the 
Bardeen-Stephen prediction, even in the limit of small ve- 
locities v < |*4f where the quasiclassical description is no 
longer valicO. Contrary, it was claimed in Rcf. [l^ that 
the additional correlation between levels (the 'two-comb' 
structure) may lead to an anomalously high vortex vis- 
cosity in the moderately clean regime in the low- velocity 
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limit. 

In our paper, we reconsider the problem of the vor- 
tex viscosity in a two-dimensional s-wave superconduc- 
tor, taking full account of the discreteness of the vortex 
spectrum and of the microscopic structure of the quasi- 
particle levels. We assume the moderately clean limit 
with an appropriate number of impurities in the core, 
such that the spectrum possesses the randomly shifted 
two-comb structure (the precise condition is given in 
Section QI). We find that in spite of the level correla- 
tions derived in Refs. p]Jl2| the vortex viscosity <j 
differ from the well-known quasiclassical resultc 
for three-dimensional superconductors. Within the same 
microscopic model we consider separately the discrete- 
spectrum and continuum-spectrum regimes of dissipa- 
tion: in both limits we arrive at the Bardeen-Stephen 
result for the dissipation, 



solutions to the Bogoliubov-deGennes equations 
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where n is the electron density, Hup ~ A 2 /Ef is the spac- 
ing between the levels in the coreta, and B is the magnetic 
field. The two-comb level structure found in Refs. [n] . |l2| 
may be classified as the ciccular unitary random-matrix 
ensemble of dimension twoE9. Therefore our findings may. 
be considered as a generalization of Wilkinson's resultsEa 
to circular ensembles. We also find the scattering time 
r for quasiparticles in the core in the s-wave approxi- 
mation: for a weak impurity strength, the effective 1/r 
is larger than the bulk scattering rate l/r„ by the log- 
arithm of the impurity strength. A similar logarithmic 
correction was previously derived in Refs. |]j7| 

The paper is organized as follows. In Section || we 
prepare for the calculation by deriving the microscopic 
Hamiltonian projected onto the relevant subgap states 
in the vortex core. In Section [II we review the results of 



Refs. |ll|,|12| on the circular random-matrix ensemble ap- 
pearing in a disordered vortex at the intermediate level 
disorder. We then describe the two limiting regimes of 
dissipation (Sect. IV), the discrete-spectrum and the 
continuum-spectrum regimes, and set up the stage for 
the calculations. In Section [v| we treat the case of the 
discrete energy spectrum (low velocities and no inelas- 
tic level broadening), while Section VI is devoted to the 



opposite continuum-spectrum limit where the levels are 
broader than the interlevel distance. Finally, we discuss 
our findings in Section VII . The sensitivity of the en- 



ergy levels to the vortex displacement is calculated in an 
Appendix. 



H A(r) 
A*(r) -H 



= E 



(2) 



where Ho = p 2 /2m — Ep. We assume an axially sym- 
metric vortex with the order parameter A(r) = A(r)e , 
where the modulus of the order parameter depends only 
on the radial component r and the phase winds with the 
angular coordinate 9. We neglect the magnetic field in 
the vortex core, assuming a large magnetic penetration 
depth A > where £ = vp/nA is the superconducting 
coherence length (here and below we choose units with 
h = 1). In the quasiclassical limit /cf£ 1, thg spec- 
trum and the eigenfunctions may be easily founda. The 
eigenvalues form a spectrum of equidistant levels 



E^ = /duo, 



(3) 



where the angular momentum [i takes half-integer values, 
fi = n + 1/2, and 



w = 



r°° Mr) -2K{r) d 
J O k F r K ^_ 

r °° e~ 2K ^ dr : 



K(r)= r^l d r>; 
Jo v F 



(4) 



(5) 



r is the distance from the core center. The basic elec- 
tronic energy scale in the vortex core takes the value 
loq ~ A 2 /ep, up to a possible logarithmic prefactor 
log(A/T) due to the shrinkage of the vorlex core at small 
temperatures T (Kramer-Pesch effectEj). The eigen- 
states ^ = (iifj,, v^) T take the form 



«V(r) = A 



■7 M -i/2(A*r)e-^-V2)^ 



(6) 



where A 2 = [ikp 1 J °° e 



-2K(r) 



dr}~ 



kp/t; is the nor- 



malization factor, and J n {x) are the Bessel functions. 

In the following, we will be interested in processes 
at energies far below the superconducting gap and thus 
project all the operators onto the subgap states ([|). It 
will be convenient to take the Fourier transform of these 
eigenvectors in the variable /i and introduce a new angu- 
lar variable <j> labelling the direction of the quasiclassical 
motion of the quasiparticle, 



J<t>/2 
-i<(>/2 



3 -- ?f ( r ')p ik « 1 



, (7) 



II. MICROSCOPIC HAMILTONIAN OF THE 
MOVING VORTEX 

Before discussing the moving vortex, let us review the 
excitation spectrum of a clean two-dimensional vortex at 
rest. The corresponding wave functions are given by the 



where we introduced the vector = kp (sin (f>, — cos <f) 
pointing perpendicular to the direction specified by the 
angle tf>, with absolute value kp. The plane- wave expo- 
nent in (^) is then k^r = kpr sin(0 — 8). This basis 
of wave functions has a very simple structure: in ad- 
dition to the phase winding exp(±i</>/2) (providing the 
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antiperiodic boundary conditions in 0), these wave func- 
tions are plane waves in the direction of the wave vec- 
tor restricted to a region of size £ around the vortex 
center. Note that ^(r) is not an eigenfunction of the 
Hamiltonian (En. Prepared in such a state at t = 0, 
the wave function will rotate in the <^-basis according to 



*(r,t) = * 



(fj — UjQt 



(r). 



We now turn to the problem of the moving vortex with 
impurities. It can be described by the time-dependent 
Bogoliubov-deGennes equations 



(8) 



where H(t) is the Hamiltonian of the vortex at the posi- 
tion r = vt with v the vortex velocity, 



m _(H + U(r) A(r-vt) 
«W - I A*(r- vt) -H -U{r) 



(9) 



and U (r) is the potential set up by the impurities. 

In the clean limit (At 1), the admixture of bulk 
states (with energies greater than A) to the vortex states 
may be neglected. We project the Hamiltonian onto 
the subgap states of a clean vortex (fj]) by substituting 



*(r,t) 



2tt 



a(0,t)^(r-vi). (10) 



Defined in this way, the amplitude a(<j>, t) has antiperiodic 
boundary conditions in 0, a(</> + 2n) = —a(<fi). The time 
evolution of the coefficients a(4>, t) obeys the Schrodinger- 
type equation 

9 ■ 8 j. f^' u ,,, [<ty T ,,,, 

t-Q^a = -iu —a + J — T < f t( ya(cp) + J — L^uty ), 

(11) 

where the kernel 



T W = (* (r- vt)| (-i^j |^(r-vt)) 

= iv<\EgV|\£y> (12) 

is produced by the vortex motion, and the kernel 

= (* (r - vt)| U{v)t z I ^(r - vt)) (13) 

is due to the impurities. 

The T-kernel may be easily computed from the explicit 
form (pp of In the limit kpt; ^> 1, the matrix element 
( |l2l ) takes the form 



b') k v 



(14) 



and we identify this term with the "Doppler shift" . 

The impurity potential U (r) is taken as a sum over 
point-like impurities, 



(15) 



Then the scattering kernel may be expressed aal 



n < 

Laa' = 2iA sin ■ 



- V i/. e - 2A: (i r *- v *i) e i (v- k «)( r *- vt ) 

(16) 

Summarizing, we end up with the equation of motion 
( |i"l| ) containing three terms: The first term describes 
the circ ular ■ ( "chiral" ) motion of the quasiparticle in the 
vortexBlLilEl The second term is the Doppler shift due 
to the vortex motion. And the third term describes the 
scattering off impurities (taking point-like impurities is 
equivalent to including only s-wave scattering). We are 
interested in the energy pumping in this time-dependent 
model. Note that the physical energy is given by the first 
and third terms in the evolution operator ( |Tl| ) but does 
not include the second term (T-term) which arises from 
the time-derivative of the basis wave functions. This dis- 
crepancy between the energy and the evolution operator 
may be resolved by an appropriate time- dependent gauge 
transformation 



a(4>) = a(4>)e' 



(17) 



This gauge transformation has a dual effect on the equa- 
tion of motion ( |TT| ) : Firstly, it makes the energy operator 
coincide with the evolution operator. The old T-kernel 
is now replaced by a time-dependent one, 



= 2w5((f> ~ 4>') w (n z x k )vt, 



(18) 



where n z is the unit vector perpendicular to the plane 
(this gauge transformation resembles the one in electro- 
dynamics replacing a static electric field with a magnetic 
field linearly growing in time). Secondly, the L-kernel is 
transformed as well and the new kernel takes the form 



Lthd,' — 2iA 2 sin • 



-2K(\ Ti - v t\) i(k^-k*)r< 



_ (19) 
which differs from Eq. (^6|) by the cancellation of the 
velocity term vt in the last exponent. 

The physical content of the two gauges may be under- 
stood in the following way: The basis wave functions (0) 
are quasiclassical plane waves (at the wave vector Uf) 
cut off by the long-wavelength envelope exp(— 2K(r)) of 
the size of order £. In the original gauge [with variables 
cl(4>, t)\, this basis was chosen by simply translating the 
basis (0) together with the vortex. In the new gauge [with 
variables a((j),t)], only the long-wavelength envelope is 
translated, without shifting the phases of the quasiclas- 
sical plane waves. 

This difference in matching phases of plane waves at 
different vortex positions produces two different descrip- 
tions of the moving vortex. In the first gauge [variables 
a((j),t)], the evolution equation ( pd| ) contains impurities 
moving with respect to the vortex (fast oscillations in 
Laa 1 )- This description is similar to the approach taken 
in Refs. [I^,|ll],[l2| a static vortex subject to moving impu- 
rities. However, in those references, the T-term was omit- 
ted. We will show below that omitting this term does not 



change the result for the dissipation in the moderately 
clean limit, thus justifying the approach of Refs. |ic|[llp4 
In the second gauge [variables a((j>,t)], there are no os- 
cillating terms in L^p (except for the slowly varying en- 
velope K(r) whose time derivative may be neglected in 
most cases) . All the oscillations in the L-term in the first 
gauge may be removed by a single gauge transformation, 
as all impurities move with respect to the vortex with the 
same velocity and in the same direction. This fact has 
not been properly taken into account in Ref. [l^, which 
has lead to an unphysical result. We shall comment in 
more detail on the derivation of Ref. |l2| in Sec. |v|; here 
we just remark that the parallel motion of the impurities 
with respect to the vortex requires special care in the 
calculations within the first gauge, but is automatically 
taken into account in the second gauge. 



III. CIRCULAR UNITARY ENSEMBLE OF THE 
QUASIPARTICLE LEVELS IN THE 
DISORDERED VORTEX 

In this section, we review the derivation of the two- 
comb spectral statistics in the vortex from Refs. [Tl| , |l~2| 
and set up the notation for the calculations in the later 
sections. For the discussion of the spectral statistics in 
this section we take the vortex at rest (v = 0). Then the 
Hamiltonian in the equation of motion ([ll]) contains only 
two terms, the kinetic term —ioj^d / d<p and the scatter- 
ing term L^p (for a vortex at rest there is no difference 
between L^p and L^). 

At not very high impurity concentration [see Eq. ( |38| ) 
below] , the scattering kernel may be approximated 
as a sum of two local terms [this approximation is due to 
the rapidly oscillating exponent in (Q)], 

N 

L<t><t>' = ^TriuJo ^[Ji^^fV ,0i+7r — ^<Va-Bt<V,4>J ) 

i=l 

(20) 

where N is the number of impurities in the core, the 
parameters 4>i specify the angular positions of the impu- 



rities, and 5, 



01, 02 



are ^-functions smeared over the width 



H and 



2) ~ (&f£) *' 2 and antiperiodically continued in 



,02+27I_ 



1+2-7T,, 



Note that 



the regularization (|20| ) as a product of two smeared 5- 
functions is important for evaluating the scattering ma- 
trix (H) belowliiy. The effective strength of the i-th 
impurity is 

J. - %A V * e ,■ ,Q.jL p -2K(u)+2ik Fri 

oj k F ri r l 

(21) 

where <&i = mVt is the Born parameter of the impurity 
(an additional imaginary unit compared to the notation 
in Ref. |l2] is due to the antiperiodic boundary conditions 
employed) . 

The matrix elements (^0|) couple only angles with dif- 
ference close to 7r. Introducing the two-component vector 




FIG. 1: Schematic representation of the evolution of the wave 
function ip(<f>) in the circular-unitary-ensemble regime. The 
two-component wave function ip(<j>) rotates with the angular 
velocity u>o and scatters off impurities. The impurity scat- 
tering is local and is described by the unitary matrices Mi, 
according to (|23|). After rotating by half a turn, the wave 
function is projected back to the origin with the help of the 
boundary condition (pj|). 



o(0) ^ 
a((j) + n] 



e [0;tt] 



(22) 



the scattering becomes local in tp. The individual scat- 
tering events in ( JTT| ) then may be integrated separately 
and formulated in terms of a boundary conditions 



where 



Mi 



iPifa + 0) = MiV>(& - 0), 



'1HJ,| 2 2J, 
1 + IJJ2 V -2J* 1-W 2 



(23) 



(24) 



The boundary condition for ip((j)) going around the half 
circle is 



V'(O) = —ia y il)(n). 



(25) 



Thus the energy levels are determined by the eigenvalues 
of the matrix —ia y M, where 

M = M N . . . M 2 M 1: 7T > <j) N > • • • > 0x > 0. (26) 

The full "scattering matrix" —ia y M is a unitary matrix 
with the eigenvalues exp(±j7r E / ujq) , so that the energy 
levels E are solutions to the equation 



cos(7r£'/ajo) = — ~tia y M. 



(27) 



The evolution of the wave function ip{4>) is schematically 
shown in Fig. |l|. 

In the moderately clean limit with a sufficiently large 
(but not too large) number of impurities in the core, the 
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matrix M (and thus also —io~ y M) is random with a uni- 
form distribution over the SU{2) group. Such a random- 
matrix ensemble is classified pas the "circular unitary" 
ensemble with dimension twcc3. The density of states 
and the level correlations may be easily computed from 
Eq. (^)--.For example, the average density of states is 
given bytHl 



(p{E)) = — sin 2 7T— 



0J 



0J 



(28) 



The spectrum consists of two combs of equidistant levels: 
A level is characterized by the comb number s = ±1 and 
its position in the comb, E± t k = ±(£0 + 2uj()k), where 
Eq is the eigenvalue of the lowest level with positive en- 
ergy (0 < Eq < ujq). The cigenfunctions can be written 
explicitly as 



^ ±!k {4>) = e^ E ^ a M^ ± {4> = 0), 



(29) 



where ip±(4> = 0) are the two eigenvectors of —ia y M, and 
the scattering matrices M^^i arc defined as 



= M in 



. Mi. 



> 



> ■ ■ ■ > fa > <f>', (30) 



where the product is taken over all impurities between 
the angles </>' and 4>. Eq. j30| ) introduces M^p for > 
</>'; furthermore, it is convenient to define the scattering 
matrix for tfi < <\>' as M^p = M^, ^. The symmetry of 
the SU(2) matrices 



(31) 



allows us to translate easily between the eigenfunctions 
for the two series of levels via 



(32) 



In our calculations, we shall need not only the proper- 
ties of the spectrum, but also some statistical properties 
of the wave functions in the regime of the circular unitary 
ensemble. At any given point <j> and for any energy level, 
the wave function tp(<p) is a spinor pointing in a random 
direction. All directions are equally probable, allowing 
us to compute equal-point correlation functions such as 



(33a) 
(33b) 



The correlations between wave functions at different 
values of cj> may be expressed in terms of the properties 
of the scattering matrices ( 30 ) . The angular correlations 
between scattering matrices decay exponentially, 



where 

1 

LOqT 



— = 16ni. 



r dr 



\J(r)\< 



(l + |J(r)|* 



(34) 



(35) 



and ni m p is the impurity concentration. In (35) we have 
assumed (for simplicity) that all impurities have equal 
Born parameters, in which case the effective impurity 
strength Jj as given by (|2l| ) becomes a function of r< and 
is denoted J{r) in ([35]). In order to derive ( [34] ) and ( |35| ) 
one takes a small increment of (M^ ^(TgM^t) in 4> and 
averages over a single impurity scattering matrix, 



-QgiMl^cTzMw) = 2n imp n 



dr; 



x '(Ml^Mj^MiM^) - (Ml ^M^,)] . (36) 



The last averaging is performed independently over 
and Mi. In averaging, we first account for the rapidly os- 
cillating off-diagonal elements in ( p4| ) and arrive at (|34|), 

For strong impurities with the Born parameter $ = 
mVi ~ 1 (as above, we assume that all impurities are 
of equal strength), the scattering time r is of the order 
of the bulk normal-state elastic mean- free time r n . For 
weak impurities with d <C 1, the integral in Eq. ( Ba ) has 
a logarithmic dependence on the Born parameterH, 



(37) 



A similar logarithmic behavior of the scattering rate was 
also derived in Refs. ^,0. The relaxation angle 4> T in fl35| ) 
plays a central role in the context of vortex dissipation, 
as it enters the expression for the friction coefficient and 
hence also for the flux flow conductivity. In Ref . [[2] this 
quantity entered the numerator in the expression for the 
conductivity, while in our results it enters the denomi- 
nator in the Bardeen-Stephen form (fir). We comment in 
detail on our disagreement with Ref. y1| in Section ^ and 
in the Appendix. 

The physical meaning of 4> T is the correlation length of 
the scattering matrix M^^i and r is the corresponding 
scattering time. We shall see later that it is the correla- 
tion function ([34]) which determines the rate of interlevel 
transitions. In principle, one can derive a similar expo- 
nential decay "in the x and y directions" [with a z re- 
placed by a x or a y in ([34|)], in which case the correlation 

length <j)i xy ' > is twice larger than <p T "in the z direction" . 
This is due to the specific form ( p4| ) of the scattering 
matrix Mj. 

In order to realize the regime of the circular unitary 
random-matrix ensemble it is necessary that <j) T <C tt, 
which defines the lower bound for the moderately clean 
regime, ojqt <C 1. However, the region of the circular uni- 
tary ensemblfedoes not extend over the whole moderately 
clean regimeE3. The additional restriction originates from 
the breaking of the instant-scattering approximation j20| ) 
in the limit of too strong disorder. Indeed, the typical 
width of the smeared 5-functions in this equation can be 
estimated as 8<j) ~ l/^/kpt; (that corresponds to an im- 
purity at the distance of the order of £ from the vortex 
center). On the other hand, the number of impurities in 
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the core is N ~ ni mp £ 2 ~ 1/wot„i9 2 . The solution for the 
spectrum discussed above is justified as long as the S- 
functions in Eq. (^0|) do not overlap, which is equivalent 
to the condition N 5(f) *C 1. Thus, the circular unitary 
ensemble is realized only in a relatively narrow range of 
disorder strengths 



In 1/1? t„ 



(38) 



[we assume that the Born parameter ^ < 1 and make 
use of Eq. (j3jj); for $ ~ 1 one should replace In 1/0 by 1]. 
For stronger disorder (i.e., larger scattering rate l/r„) 
or smaller Born parameter ■& the instant-scattering ap- 
proximation ( po|) fails and the circu lar u nitary ensemble 
crosses over to the class C ensemblecafij. A phenomeno- 
logical approach to this crossover has been discussed in 
Ref. and the results of numerical simulations are avail- 
able in Ref. E3l 



IV. TWO REGIMES OF DISSIPATION 

The time evolution of the core states is described by 
Eq. (pi]). Since the evolution operator is nonstationary 
(i.e., explicitly time-dependent), it causes transitions be- 
tween levels with different energy. In a fermionic system, 
the rate of downward transitions is suppressed compared 
to the rate of upward transitions due to the Pauli ex- 
clusion principle, leading to an increase of the average 
energy with time. The energy pumped into the system 
will finally be transferred to the thermal bath via the in- 
teraction with phonons or other soft degrees of freedom, 
thus producing a finite dissipation. 

There are two different mechanisms of dissipation de- 
pending on whether the individual energy levels can be 
resolved or notiia. If the discrete spectrum is smeared 
into a continuous one, the energy pumping can be calcu- 
lated with the help of the standard linear-response Kubo 
formulate. This regime is naturally realized when the in- 
elastic width 7i n of an energy level exceeds the mean level 
spacing luq. On the other hand, the spectrum may turn 
out effectively continuous even at 7; n <C ujq if the time de- 
pendence of the evolution operator on the right-hand side 
of Eq. ( p!l[ ) is so fast that it destroys the instantaneous 
adiabatic spectrumllS. In this case, the frequency j v of 
perturbations due to the nonadiabaticity of the spectrum 
exceeds ujq and plays the role of the effective level width. 
In the opposite case, when 7i n ,7„ <C wo> the spectrum 
is essentially discrete and— the dissipation is due to rare 
Landau-Zener transitional!! taking place when two levels 
come very close to each other. 

In a normal metal at low temperatures, the inelas- 
tic widths due to the electron-electron and electron- 
phonon interactions are given_.by j e -e ~ T 2 /Ep and 
7e-ph ~ T 3 /0 2 D , respectively^ (®d is the Debye tem- 
perature). Since 7 e - e is of the order of ujq already at 
T ~ A, it does not contribute to the level smearing at 
lower temperatures. Furthermore, for the states localized 



in the vortex core the interaction with phonons appears 
to be-,strongly suppressed compared to the normal-state 
rate£3; one of the main reasons is that the quasiparticles 
with energy (iluq <C A are composed of nearly equal mix- 
tures of electron and hole components [i.e., (u 2 ) » (v 2 ), 
cf. Eq. (Q)] with negligible net charge. Thus, at suffi- 
ciently small temperatures T <C A, the inelastic width of 
the core states 7i n <C ujq and the regime of dissipation is 
determined solely by the vortex velocity. In Ref. [lO] the 
crossover velocity vk separating the regimes with discrete 
and continuous spectra was estimated as 



A/p F 



(39) 



and we will present a microscopic derivation of this result 
below. 

For model Hamiltonians-irom the three Wigner-Dyson 
random-matrix. ensemblesc3, dissipation was calculated 
by WilkinsontScJ in the regimes of small (v vk) and 
large (u > velocities (in the limit 7i n = 0). He finds 
that if the Hamiltonian H(X) of the system depends on 
time through X = vt then the energy dissipation is de- 
termined by the variance of the (off-diagonal, i ^ k) ma- 
trix elements C(0) = (\dHik/dX\ 2 ) /lu 2 normalized by 
the mean level spacing ujq. In the continuous-spectrum 
regime specified by ti > vjc = woC _1 ^ 2 (0), the energy 
dissipation rate as given by the Kubo formula, 



WKubo = irC(0)v 2 



(40) 



is the same for the orthogonal ((3 — 1), unitary ((3 = 2), 
and symplectic (f3 = 4) ensembles and describes viscous 
damping. In the discrete-spectrum regime, at v <C Vk, 
the result depends on the level-repulsion parameter /3, 
which determines the probability P(e) oc e 13 to find two 
levels at a distance e <C u>o- In this situation the energy 
dissipation rate is given by the expressionllj 



(41) 



hence the dissipation is superohmic for the Gaussian or- 
thogonal ensemble, while for the Gaussian unitary en- 
semble it remains ohmic with W% = 7rC(0)w 2 , exactly co- 
inciding with WKubo despite a very different mechanism 
of dissipation. 

In the next sections we calculate the dissipation rate 
for the vortex motion in th£ regime of the Koulakov- 
Larkin 'two-comb' spectrun£3, both in the limits of small 
and large velocities, thus extending Wilkinson's consid- 
erations to the case of circular unitary ensembles. 



V. DISSIPATION IN THE DISCRETE 
SPECTRUM (LANDAU-ZENER REGIME) 

In this section we calculate the dissipation for the two- 
comb level structure described in Section III at small 
velocities. We will work in the 'tilded' basis a((j>, t) intro- 
duced in Eq. (|l7]) , where the energy coincides with the ex- 
pectation value of the evolution operator (p^), with Ta,m 
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and L,i,th ' substituted by their tilded counterparts (|1^ 
and (pL9|) . In the low- velocity limit, the energy is pumped 
into the system when two levels come very close to each 
other and a nonadiabatic Landau Zcncr transition be- 
comes possible. Due to the symmetry of the two-comb 
spectrum the rate i? LZ of such transitions is the same for 
each neighboring pair of levels. We will calculate R LZ by 
considering the lowest level with positive energy e and its 
mirror image with energy — e. For simplicity we assume 
that the vortex moves in the x direction, v y = 0. 

Following the logic of Ref. [18], we diagonalize the 
Hamiltonian at t = and restrict it to the 2x2 ma- 
trix involving the pair of states considered, 



H(t) 



e + hut hi 2 t 
h\ 2 t —e — hut 



(42) 



where hij are the elements of the matrix dH/ dt (neglect- 
ing the quadratic term oc (d 2 H/dt 2 )t 2 is justified in the 
Landau-Zener regime since the duration of a nonadia- 
batic transition is proportional to the gap of the avoided 
crossing and is small at v <C vk)- The instantaneous 
adiabatic spectrum takes the form 



E(t) = ±^A 2 z + A 2 (t-t, z ) 2 , 



(43) 



where 



e 2 \hi 2 
A 2 



A 2 = h\i 



112 



ehu 



A 2 



(44) 



Equation (43) describes an avoided crossing with the 
minimal distance (Landau-Zener gap) 2A LZ between the 
spectral branches realized at t = f LZ . The probabil- 
ity of the Landau-Zener transition at such a crossing is 
exp(— ttA 2 z /A). The mean rate of transitions is given by 



5(i LZ )exp ( -tt^ 



(45) 



where the role of <5(t L z) is to count each avoided crossing 
once. The average in Eq. ( |45| ) is taken over the distri- 
bution of the parameters e, hn, and describing the 
avoided crossing. The energy e is expressed through the 
transfer matrix M according to Eq. (p7|). The coefficients 
/in and hi2 are the matrix elements of the perturba- 
tion ( |l8| ) over the exact wave functions ip±(4>) =_0±,o(0) 
[which depend on the trajectory M^o via Eq. (|29|)1, 



hn = loqUfv / — ^t(0)cr z V'+(0) cos</>, (46) 
Jo T 

hi2=uj k F v — ^\.{4>)ct z ^-{4>) cos0. (47) 



We come to the crucial point: The quantities e, hn, 
and hi2 have different dependencies on the transfer ma- 
trix -M^.o- In the moderately clean limit, when the num- 
ber of impurities in the core is sufficiently large and 



<p T <C 1, the matrix M^o performs many rotations over 
the SU{2) group. Therefore, we conclude that (i) e, hn, 
and hi2 are uncorrelated, and (ii) the distribution of hn 
and hi2 is Gaussian. We further calculate the variances 
of hn and hi%. For (hn) we obtain with the help of 

Eq. m, 



(h 2 n) = (u k F v) 2 I ^^costf<We-l*-*'l/*' 



r7r d(j) d(j>' 



(48) 



Integrating over d>—(f>' and averaging over the spinor ip((j)) 
ind 

^((^ + ) 2 ) = ^- (49) 



according to (33a) we find 
(hli) 



(uj k F vy 



3tt 



Analogously, with the help of Eq. (33b) we obtain 



(\hi 



(u; k F v) 2 



— (\ip+<J x ip 



* |2 
+ 1 



20r 

3tt ■ 



(50) 

In the same way it can be shown that (/in/112) = 0, 
thus proving the statistical independence of hn and hi 2 . 
We conclude that hn, and Qhi2 are independently 

distributed with the same distribution 



T(h) 



1 



2</> T uiokpv 



exp 



3nh 2 



2(j) T (u>okFvy 



(51) 



The distribution function for e is provided by the density 
of states (H||), in the limit e < uj q , V{e) = 2tt 2 s 2 /lu^. 

Averaging Eq. (^B| ) over e, hn, and hi2 with the dis- 
tributions V(e) and V(h), we obtain for the mean rate of 
Landau-Zener transitions 



R, Z = kJ f^- (52) 



The energy dissipation rate is given by W 
with the vortex viscosity 



uj R lz = rjv 2 



(53) 



As defined by Eq. (J53J), 77 is a two-dimensional viscosity. 
In a layered superconductor, it determines the friction 
force F = — rjv exerted on the vortex from excitation of 
quasiparticles within one layer. 

The statistics (^TJ) of the matrix elements of dH/dt 
determines the sensitivity of the spectrum to the vortex 
motion and allows one to find the critical velocity vk 
separating the regimes of discrete and continuum spectra. 
The discrete spectrum can be resolved if the change of the 
Hamiltonian during the time ujq 1 , 5H ~ h/cjQ, is smaller 
than luq. Taking for h the dispersion of the distribution 
(pl|), we find for the critical vortex velocity 



v K 



(54) 
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coinciding with the estimate (|39|) found in Ref. |1Q. 

Comparing with Ref. [l2| we find that our result (|5S 
for the viscosity is smaller by a factor (f) 2 . and the expres- 
sion (HJ) for vk is larger by the factor 1/4> T . In Ref. [IJ, 
following the analysis in Refs. it was assumed that 

the dissipation is produced by impurities moving with re- 
spect to the vortex, equivalent to neglecting the Doppler- 
shift term T^i in Eq. ([ll]) . In order to identify the origin 
of the discrepancy in the results and verify the validity 
of omitting the Doppler shift, we reconsider below the 
treatment of Ref. [l2] — we will see that the source of the 
disagreement is not in the neglect of the Doppler-shift 
term but in the neglect of cross-correlations between the 
motion of different impurities. To demonstrate this, we 
adopt the approach of Ref. [l2| (omitting the Doppler- 
shift term) and re-calculate_the correlator of the spatial 
gradients of the matrix MtH 



C 



tr 



dx J 



(55) 



from which the viscosity follows via r\ = 2C/n. Compar- 
ison between Eqs. (55) and ([l0|) shows that the quantity 
(2/ir 2 )C plays the role of the correlation function C(0). 

The matrix M specified in Eq. ( p6| ) is a product of 
transfer matrices of individual impurities. Therefore, 
there are two contributions to the correlator (^5|), one 
originating from derivatives over the coordinates of the 
same impurity (diagonal part) and a second one from 
derivatives taken on different impurities (off-diagonal 
part), 



lr^/ dM dM\ 
Cdia g = I ^(tr — tr— J 

dM\ 

dxi dxj I 



dxi dxi / ' 

lv^/ dM dM\ 
Co ff -d iag = I ^(tr — tr — ) 



The diagonal part is given b 



yB 



irk 2 F 



(56) 
(57) 

(58) 



The cross term, expressing the correlated nature of the 
motion of the impurities with respect to the vortex, has 
been missed in Ref. O. Its calculation is presented in the 
Appendix and the result takes the form 



^off— diag — 



TTfcf, 



4^.(1 + #)" 



(59) 



We thus see that the two contributions Cdiag and C ff-diag 
nearly cancel each other and the net sensitivity of the 
spectrum to the vortex motion appears to be significantly 
lower than that found in Ref. 



C = 



nk%4> r 



(60) 



where we employed the condition <j) T <C 1 of the moder- 
ately clean limit. Using Eq. (|55|) one recovers the result 



( p3| ) for the vortex viscosity. The above consideration not 
only corrects the result of Ref. in addition, it serves 
as a microscopic justification of the model adopted in 
Refs. ^p|Jll| , [l^ , where the dissipation is due to impurities 
moving through the core and the Doppler-shift term T 
in Eq. ([Til) is neglected. 



VI. DISSIPATION IN THE CONTINUUM 
SPECTRUM (KUBO FORMULA) 

If the vortex velocity v exceeds vk the discrete spec- 
trum is smeared and transitions between non-nearest lev- 
els become possible. In this limit, the energy dissipation 
can be calculated with the help of the Kubo formula; as 
m Sec. we will use the 'tilded' basis (|l7j). Neglecting 
the slow time dependence of L^p , we write the evolution 
operator as H(t) = Hq + H, where Ho is the Hamiltonian 
of the vortex at rest and I — T/t takes the form 

I u > = 2irS(cj)- <j>')I(4>), I{4>) = u k F vcos(j), (61) 

where we again assumed that v y = 0. The energy dissi- 
pation rate 



(62) 



is calculated as the linear response to the term It with 
the help of the Kubo formula 



W = i 



t'dt'(0\[l{t'),l(0)}\0), 



(63) 



where 



T(i)= / X /(#t^ t )f ( ^t)^ (64) 

and ^f((f>, t) is a second-quantized operator in the interac- 
tion representation. The latter can be rewritten in terms 
of the (exact disorder dependent) eigenfunctions a n (<p) 
and eigenvalues E n of the Hamiltonian Hq, 



-iE n t 



(65) 



with c n the corresponding annihilation operators. Aver- 
aging over the initial state of the vortex at rest (so that 
a n {4>) coincides with a n (0)) one arrives at the expression 



(27T) 

x ak{(t>)al{(j)')a m {(j)')a* m {(j))e 



i(E h -E m ){t'-t) 



, (66) 



where = (c^Cfe) is the distribution function at the 
energy Ek- This expression can readily be represented 
(cf., e.g., Ref. 151) in terms of the Green functions 



G 



R,A 1 



E-E n ±iS' 



(67) 
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W =2 



(2n) 



dE dn{E) 
2^ dE 



mm 



(68) 



The above formula for the dissipation rate is not re- 
stricted to the regime of the circular unitary ensemble 
but is valid in the whole moderately clean limit regime 
as long as v 3> vk- In particular, one can easily recover 
the result of the standard quasiclassical analysis by as- 
suming the T-approximation: within this approach the 
Green function is diagonal in the momentum representa- 
tion, 



<g*' a (m)> = 



1 



E - Lu of i ± j- 



(69) 



and evaluating the integral in (|6^) one finds the energy 
dissipation rate W 
efficient 



rj T v with the vortex viscosity co- 



-k 2 F uJoi 



(70) 



The T-approximation is the simplest approach to the 
problem where all spectral correlations are neglected. Be- 
low, we will calculate the vortex viscosity for the case of 
the Koulakov-Larkin two-comb statistics and check how 
the presence of the two-comb correlations in the energy 
levels modifies the result ([70|). To this end we rewrite 
Eq. (|68| ) in terms of the wave functions tp s k (<fi) introduced 
in (|29^ 

s,s' ,k,k' 

x nfflW)^, ttMM •V4(</>'K</W(0')> 1 

where $ = E(t — t') — E s kt + E s ik>t'. We first re-express 
the wave functions ip s k(<t>) via (|2^) and perform the sum- 
mation over k and k' according to the summation formula 



^ e ^ s ,ot Y S(t-rmr/uo) (72) 



ra— — oo 



(which produces a double infinite sum of ^-functions of t 
and t'). We cut off the integrals over t and t' at some time 
scale smaller than ir/uo — such a cut-off is equivalent to 
assuming a smearing of the energy levels with widths 
larger than ujq. Under this assumption, only one term of 
the double infinite sum survives with t = if = {4>— <t>')/u>o- 
Next, the integrations over t, t', and E can be trivially 
performed. Finally, we sum over s and s' and arrive at 

(^(0) Ml ^M^o ^(0)^1(0) M^M^o Vv(0)> 

= ti{M^a t M^a t ). (73) 



This average is calculated with the help of Eq. ([34]) yield- 
ing 



dcj)d(j)' cos0 cos^' tr(M^, ,a z M^' t ^a 



-kpUjQT. (74) 



k'p 

47T 



The form of the above result coincides with ( |70[) cal- 
culated within the r-approximation. The microscopic 
expression for the elastic relaxation time r is given by 
Eq. (|3"5|). In case the inelastic relaxation time is shorter 
than r, it substitutes the latter in Eq. (|74|). 

To conclude this section we mention that the same re- 
sult (fr3) can be obtained in the model where the Doppler- 
shift term in Eq. (|ll]) is neglected and the dissipa- 
tion is due to the motion of impurities with respect to the 
vortex. Within that model the energy dissipation rate is 
calculated as a linear response to a change in the impurity 
positions. The derivation formally repeats the one pre- 
sented above but with the operator I^i — v dL^ /dx. 
Performing the same manipulations leading to Eqs. j7l| ) 
and (frl) we obtain 



V 



1 

47T 



tr- 



dW dM 
dx dx 



(75) 



where M is the product of scattering matrices defined 
in (p6j). The above correlation function for the random 
SU(2) matrix M can be related to the correlation func- 
tion ( J55| ) and is equal to 8C Therefore, Eq . (|75|) exactly 
coincides with its low- velocity analogue (55) and repro- 
duces the result (frl). 



VII. DISCUSSION 

The main conclusion of this paper is that the Bardeen- 
Stephen expression for the flux flow conductivity is ex- 
tremely insensitive to the details of the level correlations 
in the vortex core. We calculated the energy dissipa- 
tion rate for the case of the Koulakov-Larkin two-comb 
level structure when the spectral correlations are most 
pronounced. Such a statistics, classified as the circu- 
lar unitary ensemble of dimension two, is obtained for 
layered superconductors within the region of the moder- 
ately clean limit specified by Eq. (Q). We found that 
the vortex viscosity is the same for the regimes of dis- 
crete (v < vk) and continuum (v > Vk) spectra and 
coincides with the result (^o|) of the phenomenological 
r-approximation. The viscosity i] determines the flux- 
flow dissipative conductivity via the standard relation 
Oxx = ecq/irB. Assuming a cylindrical Fermi surface, we 
arrive at the result (|l|) Pjj4iWWri previously via several 
quasiclassical approachcsaBal§uQ. 

Despite we found the same result for a xx in the limits 
of small and large vortex velocities, it is worth emphasiz- 
ing that the physics of energy dissipation is quite differ- 
ent in the two limits. In the Kubo regime [v > Vk) the 
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energy is pumped continuously, whereas in the Landau- 
Zener regime (v < vk) it is absorbed in discrete portions 
of size u>q. The equivalence between the r-approximation 
and the result of the exact microscopic treatment is not 
very surprising in the Kubo regime, as in this limit it 
is determined by the net sensitivity of the spectrum to 
the vortex displacement rather than by interlevel corre- 
lations. On the other hand, in the Landau-Zener regime, 
this equivalence is a matter of coincidence: it relies on 
the fact that the vortex Hamiltonian belongs to the cir- 
cular unitary universality class characterized by a level 
repulsion parameter (3 = 2. 

The other outcome of the present work is that it pro- 
vides a justification for neglecting the Doppler-shift T- 
term, E q. that was implicitly assumed in earlier 

papersE2ltlHld Without the T-term, the system is equiv- 
alent to a vortex at rest with impurities moving through 
its core. The dissipation is then related to the change of 
the level positions due to the motion of the impurities. 
On the other hand, the T-term describes the electric field 
in the core, which is the source of energy dissipation in 
the Bardeen-Stephen model — omitting this term is a bit 
confusing. Nevertheless, our analysis indicates that the 
vortex viscosities calculated with and without this term 
coincide in the moderately clean limit. One may expect, 
however, that a correct treatment of this term is cru- 
cial in the calculation of the Hall conductivity and in the 
superclean limit. 

In this paper we considered a strictly two-dimensional 
superconductor (or layered superconductor with negligi- 
ble interlayer coupling). In the case of the Koulakov- 
Larkin circular unitary ensemble statistics, the effect of 
interlayer hopping is small provided that the effective 



quasiparticle temperature in the core, T corc , is smaller 
than u)o(m c /m a b), where m c /m a ;, is the effective mass 
anisotropy. This condition ensures that the interlayer 
hopping amplitude is always smaller than ujq and can be 
neglected. This behavior is to be contrasted with the 
one for the class C random-matrix statistics realized in 
the moderately clean case for weak impurities (i9 -C 1) 
and in the dirty (At <C 1) limit. In that case, the inter- 
layer coupling leads to tunneling from the n-th level in 
one layer to the m-th level (with m =/= n) in the adjacent 
layer, thereby opening an interlayer channel for energy 
dissipation. The competition between the intralayer and 
interlayer channels may lead to a strong deviation from 
the Bardeen-Stephen formula ([[]) aneLeven to hysteretic 
behavior of the current- voltage curveE3. On the contrary, 
the two-comb spectrum is absolutely rigid: the n-th level 
in one layer can have an avoided crossing only with the 
n-th level in the adjacent layer and the tunneling between 
these two levels does not result in dissipation. 
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APPENDIX 

Here, we calculate the off-diagonal correlation function C ff_diag defined in Eq. (p7|). To this end we divide the 
angle interval [0, 7r] into many small pieces {cf>( k ~ 1 ' > , 4>^>] of width Stp — > so that each piece contains one impurity at 
maximum. The transfer matrix of the fc-th interval is hence either M^k) M-i) = 1 (no impurities) or M^k) ^(fc-i) = Mj 
(if the angle of the i-th impurity fa £ [i^* -1 ), </>(*)]). Then 

,-, -v ; dM^k) t<j) (k-i) 

f-off-diag — J 2.^1 \ trM 7r,<M fc) ~q^. M 4>( k - 1 )fl x i rM 7r,0(p) -Q- M 0(p-i>,O 



k^p 

where R = Mm) >0 M v an d 



kjtp 

lj2( txRY ( k ) x tTRM l^^ f (p) M <t>^ ( 76 ) 



Y(k) = — 4g — mJw,^-!,, yip) = «Jw,^-d— *^ — ■ ( 77 ) 

The representation ( |76|) is suitable for averaging over disorder, since one can independently average the matrices R, 
Y, Y, and M^p-i) ^k) . The statistical independence of Y(k), Y(p), and M^p-ii^w follows from the fact that the 
intervals [0( fc-1 \ </>(% (p ~ 1} ], and [^p- 1 ) , do not overlap. As we will see below, the correlator (|7^) is 

essentially nonzero at \<p^ — <f>^ \ ~ <p T . The matrix R, which couples (j)^ and (j)^ through the point <p = n, is the 
product of a large number of matrices Mj and looses all correlations within the interval [</>( p ~ 1 - ) , <f>^]. 
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In calculating (Y(k)) and (Y(j>)) with the help of Eqs. (B4T) and (|2l|) only the fast phase of Ji should be differentiated, 



(Y(k)) = -(Y(k)) =ia z 2k F cas<j) {k) 66 I 2nrdr „ ,„\. „ = i<r z k F ^ cos 6<< k \ (78) 

JO (1+KWD 9r 

with the last relation following from the definition ( p5| ) of the angle coherence length 6 T -C 1. Taking the continuum 
limit 56 — > we obtain 



2 



r — F 

Ooff-diag - 



j2 I dfa d6 2 cos 6i cos <fo (tr i?cr z x tr RM\ ,a z M^.^) . (79) 
Averaging over M^ lt ^ 2 according to Eq. (|34j) and integrating over (</>i + 6 2 )/2 we find 

Coff-diag = 7-f (tr Ra z x tr i?a z ) / d6 cos 0, (80) 



where the upper limit is substituted by infinity due to the fast convergence of the integral. Finally, averaging over R 
uniformly distributed over the SU{2) group and evaluating the remaining integral, we arrive at Eq. (|59|). 
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